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Abstract 

We present results from a Monte Carlo simulation of the Nambu - Jona- 
Lasinio model, with continuous SU(2) SU(2) chiral symmetry, in four Eu- 
clidean dimensions. Different model equations of state, corresponding to dif- 
ferent theoretical scenarios, are tested against the order parameter data. The 
results are sensitive to necessary assumptions about the shape and extent of 
the scaling region. Our best fits favour a trivial scenario in which the loga- 
rithmic corrections are qualitatively similar to those predicted by the large 
Nf approximation. This is supported by a separate analysis of finite volume 

corrections for data taken directly in the chiral limit. 
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I. INTRODUCTION 



The Nambu - Jona-Lasinio (NJL) model [Q is a quantum field theory of fermions in which 
the only interaction is an attractive contact interaction between fermion and anti-fermion. 
Its Lagrangian density is 
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1^ = i'api.f? + rn)^ap - -^TT- {i^api^apf - H'apl^Tpqi^aqT , (1-1) 

where Tpq are Pauli spin matrices with indices p, q running over an SU(2) internal symmetry, 
and a runs over Nj distinct fermion flavors. For bare mass m = there is an SU(2) ® SU(2) 
chiral symmetry 

ijL^UijL ; ^R^V^R, (1.2) 

where U and V are independent global SU(2) rotations. The original motivation for studying 
the model ( |1.1|) is that for coupling g"^ larger than some critical it gives a qualitatively 
good description of dynamical chiral symmetry breaking in strong interaction physics, with 
order parameter S oc (ipip)- At the large coupling strengths required for chiral symmetry 
breaking, the fermions interact by exchange of composite scalars and pseudoscalars. In 
the chiral limit the pseudoscalars are Goldstone modes, and are associated with physical 
pions. More recently, NJL models have been proposed as a description of the Higgs sector 
of the standard model 0, in which the Higgs scalar appears as a fermion - anti-fermion 
bound state, and the Goldstone modes are absorbed into longitudinal gauge boson degrees 
of freedom. 

In weak coupling perturbation theory (WCPT), the NJL model is non-renormalisable 
for d > 2, being plagued by quadratic divergences in four dimensions. However, a pertur- 
bative expansion made about the critical point g^ = g1 in powers of is known to be 

renormalisable for 2 < d < A P| 0], and in c? = 4 suffers from only logarithmic divergences, 
implying that a moderate separation between UV and IR scales can be made. Because the 
UV scale can never be completely eliminated from the theory by renormalisation without 



2 



rendering the interaction strength vanishing, the model is trivial for d = 4; it shares this 
property with the model more conventionally used to describe the Higgs sector, namely the 
0(4) linear sigma model, which contains elementary scalar fields. Both models are believed 
trivial, and hence provide effective field theories only for scales <^ A, the UV scale. At scales 
~ A, logarithmic scaling violations will manifest themselves, which in turn necessitate the 
introduction of higher dimensional operators corresponding to new physical information. 
Taking these new operators into account, then in the sense that the number of undeter- 
mined parameters is the same for both NJL and sigma models, the two models have equal 
predictive power for the standard model 

If we formulate definite lattice (ie. bare) actions based on these models, however, an 
approach perhaps more familiar in condensed matter physics, then triviality may manifest 
itself in distinct ways. For both models the upper critical dimension is four, ie. the symmetry 
breaking transition is described by mean-field (Landau-Ginzburg) critical exponents with 
logarithmic corrections. In the sigma model, deviations from mean field behaviour are well- 
described by WCPT, since the renormalised coupling is bounded by the fixed point value 
Qfi = 0{4 — d) and is small as d ^ 4. Using the renormalisation group ||^ 0|, the following 
prediction for the equation of state, relating the order parameter S to the reduced coupling 
t = — and the symmetry-breaking field m (we retain the notation of the NJL model) 
is derived: 

Note that E is measured in units of the cutoff, so that the argument of the logarithms diverges 
in the continuum limit. By contrast, symmetry breaking in the NJL model occurs for large 
coupling strength, and WCPT is inapplicable. An alternative expansion parameter, 1/Nf, 
can be used ||^ ; at leading order this yields a qualitatively different form for the equation of 
state 

m = AtE + 5E^ In (^^ j . (1.4) 



The different behaviour can be traced to the fact that in the sigma model the scalar excita- 
tions are elementary, and the associated wavefunction renormalisation constant Z perturba- 
tively close to one, whereas in the NJL model the scalars are composite, so that Z vanishes 
in the continuum limit for < 4 where the model is renormalisable, and Z oc l/ln(^) for 



(i = 4 [H hence we expect ( |1.4| ) to hold beyond leading order in 1/A^/. 

The form (|1.3| ) in the sigma and related scalar models has been subject to extensive 
analysis via Monte Carlo simulation (two recent examples are refs. |]10| [0), which ap- 
pear to provide support, although logarithmic corrections are notoriously hard to isolate 
numerically. For the fermionic case, a simulation of the d = 4 Gross-Neveu model, with 



discrete chiral symmetry, has shown support for the form (|1.4|) ||T2| . The issue of which is 
the appropriate equation of state is important for several reasons; for instance, it informs 
the study of the inherently non-perturbative chiral symmetry breaking transition in lattice 



QED4. Some studies |jT3| have assumed a scenario of triviality for this model based on 



a form similar to ( |1.3| ), which may not be appropriate. Another possibility is that different 
lattice implementations of the chirally symmetric four-fermi interaction, namely the "dual 
site" formulation ||T5| used in |T^ and in this study, and the "link" formulation used in a 
Monte Carlo study of the U(l) NJL model in [|r^, may lie in different universality classes 
with distinct patterns of logarithmic scaling violations. Studies of the corresponding lattice 
models in three spacetime dimensions ("dual site" in ||^ and "link" in Ijl^); in which the 



two implementations yield radically different behaviour, support this scenario. 

In this paper we extend the work of [l^ with a numerical simulation of the SU(2)(8)SU(2) 



NJL model, for a small number of flavors A^j = 2. The main motivation for this choice of four- 
fermi model, apart from the phenomenological reasons cited above, is that it has a maximal 
number of Goldstone excitations, all of which contribute to the 0{1/Nf) quantum corrections 
to the equation of state. We hope to discriminate between the two forms ( p..3| ) and (p.. 41), 



and in particular to see if ( |1.4| ) persists beyond leading order in 1/Nf. In Sec. II we outline 



the lattice formulation of (|1 . 1|) using staggered fermions, and discuss its interpretation in 



terms of continuum four-component spinors. It turns out that in our simulation each lattice 



species describes eight continuum flavors, and hence the Nf = 2 model requires 



lattice flavors. This necessitates the use of a hybrid molecular dynamics algorithm |T8 
||19|| . Extensive studies on small systems reveal the systematic error in the algorithm to 
be O(A^^Ar^), where At is the timestep used in the evolution of the hybrid equations of 
motion; hence the systematic error can be kept under control. Of course, we are still left 
with an unquantified error due to the use of a fractional A^, which implies a non-local action, 
but have no reason to suspect there is a problem. 

In Sec. Ill we confront data taken on a 16^ lattice with a range of bare mass values 
in the critical region 1/(7^ ^ 0.5, with five prototype equations of state refiecting various 
theoretical prejudices, including the simplest mean field ansatz, and a non-trivial contin- 
uum limit described by non-mean field critical indices, as well as ( |1.3| ) and ( p..4|) . No form 



satisfactorily fits all the data; this is because there is a finite scaling window where the 
correlation length is sufficiently large for a continuum description to apply. We show that 
the most successful equation of state ansatz, in the sense of small depends on which 
assumptions are made about the extent of the scaling window: if we include all mass values 
from a narrow region about the critical coupling then forms similar to ( |1.3| ) are the most 
successful, whereas if only small mass values from a wide range of couplings are used, then 
( |1.4|) is the most appropriate form. 

In Sec. IV we examine data taken at vanishing bare mass, which is possible in four-fermi 
models formulated via the dual-site approach. Due to the masslessness of the pion in this 
limit the data is subject to a significant finite volume effect, which we attempt to account 
for using a formula first derived for the sigma model We find that the variation of 



the finite volume correction with can best be explained by assuming an equation of state 
( |1.4| ), together with a wavefunction renormalisation constant Z cc 1/ln in other words, 
it is the finite volume corrections that seem most sensitive to the composite nature of the 
scalar excitations. The distinct triviality scenarios in NJL and sigma models ||^ [|T^ are 
rephrased in terms of the pion decay constant /^r. In Sec. V we present brief conclusions 
and directions for further work. 



II. LATTICE FORMULATION OF THE MODEL 



The lattice action we have chosen to study is written 

N 2 AT 

s^os = T.T.^ii^)iM^M)-i^M + — E (^'(5^) + m-m) , (2.1) 

where are bosonic pseudofermion fields defined in the fundamental representation of 
SU(2) on lattice sites x, the scalar field a and pseudoscalar triplet vf are defined on the dual 
lattice sites x, and the index a runs over species of pseudofermion; the relation between 
and the number of continuum fiavors Nf will be made clear presently. The fermion kinetic 
matrix M is the usual one for Gross- Neveu models with staggered lattice fermions [Q, 
amended to incorporate an SU(2) ® SU(2) chiral symmetry: 



(2.2) 



where m is the bare fermion mass, the SU(2) indices p, q are shown explicitly, e{x) = 
^_iyi+x2+x3+x4,^ and the symbol {x,x) denotes the set of 16 dual sites x adjacent to the 
direct lattice site x. The 3 Pauli matrices r* are normalised such that tr(rV-') = 25*-'. 
Integration over the pseudofermions yields the following path integral: 

2N 

In terms of Grassmann fermion fields X? C) C ^ can be derived from the equivalent action 



Z (X j DaDn det^ (M^M) exp ^(a^ + vf.Tr) j . (2.3) 



a 9 X 

The auxiliary fields a and vr can then be integrated out to yield an action written completely 
in terms of fermion fields: 

^fer = EE^"(^)(^ + 

xyXa (y) 

a xy 

2 
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E 



1 . . .A A / 1 



8A^ 



E E Xa{x)Xa{x) - E Ye ^ ^ Xa{x)€{x)fXa{x) 
" {x,x) I i=\ \ a 



(2.5) 



with ^ xy shorthand for the free kinetic term for staggered lattice fermions. So, we see 
that the action ( |2.1| ) describes 2N flavors of lattice fermion, but with a flavor symmetry 
group U(A^) ^\J{N). The interaction term in (|2.4|j2.5| ) has been normalised so that the gap 
equation coming from the leading order approximation agrees with that derived for 
lattice Gross-Neveu models having Z2 Q or U(l) ||2D| chiral symmetries. 



Next we identify the global SU(2)®SU(2) chiral symmetry of the lattice model. This is 
most manifest in the form (|2.4| ). Let Ve{x),Vo{x) be even and odd site projectors defined 
by 

n/„(x) = ^(l±e(x)). (2.6) 

Then, noting that Ve^Vo = Ve^ etc, we find that (|2.4|) is invariant, in the chiral limit 
m — >■ 0, under the combined transformation: 



$ = (al + m.r) f-^ (2.7) 

where U, V are independent SU(2) rotations. Note that (^^ ) is a symmetry because $ is 
proportional to a SU(2) matrix, and the auxiliary potential proportional to tr$^$. This 
property does not generalise to larger unitary groups. The symmetry (|2.7|) is broken explic- 
itly by a bare fermion mass, and spontaneously by the condensates (xx); (CO- 

Now let us consider the continuum fiavor interpretation. It is well known that four- 
dimensional staggered lattice fermions can be interpreted in terms of four fiavors of Dirac 
fermion, which decouple at tree level in the long- wavelength limit. One way of seeing this 



is via a transformation to fields q, q defined on the sites y oi a lattice of spacing 2a [21 



Q^'^iy) = ^ E rTx(^; y) ; = ^ E x{A; y^f, (2.8) 

where a, f] = 1, . . . , 4 run over spin degrees of freedom, a,b = 1, . . . , 4 over fiavor, and A 
is a four-vector with entries either or 1 which identifies the corners of the elementary 



hypercube associated with site each site x on the original lattice is thus mapped to a 
unique combination [A]y). The 4x4 matrix Va = li^l2'^lt''^lt^ ■, where 7^ are Dirac 
matrices. The action written in the form ( p.4[ ) contains interaction terms of the form axXi 
TT-x^x etc. In the g-basis the a and vr fields are not all equivalent. If we label the dual site 
(xi + |, . . . , 0:4 + |) by {A; y), then for A = (0, 0, 0, 0) the interaction reads 

S^nt{0■, y) = a(0; y)q{y){^ ® ^Hv) + ^^(0; y).g(2/)r(75 ® iMv)^ (2-9) 

where the first matrix in the direct product acts on spin degrees of freedom, and the second 
on flavor. Equation ( p. 91 ) resembles the continuum interaction between fermions and the 
auxiliary scalars, except that the pseudoscalar current is non-singlet (but diagonalisable) in 
flavor space. However, for = 1, A^^-^fj^ = 0, say, the interaction is more complicated: 



+i7T{fi;y). 



q{y)q{y) + (q{y)q{y) + q{y){l^,l^ ® 75*7^ 

q{y)r{i5 ® iDqiy) + (g(i/)r(75 ® -fl)q{y) + q{y)T{^f, ® %)(i{y) 



(2.10) 



where is the forward difference operator on the y-lattice. Hence in addition to continuum- 
like interactions there are extra momentum-dependent terms coupling the auxiliary flelds 
to bilinears which do not respect Lorentz or flavor covariance. For A vectors with more 
non-zero entries there are still more complicated interactions containing more derivative 
terms. 

One might naively claim that the terms in A^a, A^tt are 0(a) and hence irrelevant in the 
continuum limit. This ignores the fact that a and tt are auxiliary flelds which have no kinetic 
terms to suppress high momentum modes in the functional integral. The non-covariant terms 
will thus survive integration over a, vr to manifest themselves in ( p.5|) . A similar phenomenon 
happens in lattice formulations of the NJL model in which the interaction term is localised 
on a single link |T^ [|22[ . The effects of the non-covariant terms in fact depend on the details 
of the model's dynamics. If long-range correlations develop among the auxiliary flelds, ie. 
if an effective kinetic term is generated by radiative corrections, then the term may 
become irrelevant - this appears to be the case in two |]15| |^ and three dimensional 
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models, where in each case there is a renormahsable expansion available. In four dimensions 
the NJL model is trivial: hence no continuum limit exists. Even in this case, though, we 
expect correlations to develop between the auxiliary fields in the vicinity of the continuous 
chiral symmetry breaking transition. In this case the unwanted interactions would manifest 
themselves as scaling violations if an interacting effective theory is to be described |^. 

From the arguments following ( p.5|) and (|2.8| ) we now state the relation between the 



parameter and the number of continuum flavors Nj: 

Nf = 8N. (2.11) 

The extra factor of 2 over the usual relation for four-dimensional gauge theories arises from 
the impossibility of even/odd partitioning in the dual site approach to lattice four-fermi 
models: in other words the matrix M contains diagonal terms which are not multiples 
of the unit matrix. In this paper we wish to study a small number of continuum flavors 
Nf = 2, which necessitates a fractional = 0.25. This can be achieved using a hybrid 
algorithm, which produces a sequence of conflgurations weighted according to the action 
( [^.3D by Hamiltonian evolution in a flctitious time r |T^, the flelds' conjugate momenta 
being stochastically refreshed at intervals. The advantage of this approach is that in the 
form ( p.3|) the variable can readily be set to a non-integer value (though there is no longer 
a transformation to a local action of the form (|2.5| )); the cost is that the simulation must be 
run with a discrete timestep At. In principle several values of At must be explored, then 
the limit At taken. We have chosen to implement the "R-algorithm" of Gottlieb et al 
T9| , for which the systematic error is claimed to be O(Ar^). 

We have tested the R-algorithm by extensive runs on a 6^ lattice with parameters 1/ g"^ = 
0.56 and mass m = 0.01, the coupling being chosen to lie well into the broken symmetry 
phase according to the leading order gap equation (see below). We tested models with 
= 0.25, 1, 3 and 6, with timesteps At = 0.4 (forA^ < 1), 0.2, 0.1 and 0.05. For integer A^ 
we were also able to run directly in the At = limit using an exact hybrid Monte Carlo 
algorithm H @]. We ran for either 20000 (A^ = 3,6) or 40000 (A^ = 0.25, 1) Hamiltonian 
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trajectories between momentum refreshments, the trajectory lengths being drawn from a 
Poisson distribution with mean 1.0. 

For illustration we present results for two local quantities, the expectation value of the 
scalar field (a) = S, and the energy density e given by 

YM-^ , -M-^ A. (2.12) 

X I 

The results are shown in Figs. ^ and |^, together with quadratic fits of the form 

S(Ar; Ar) = Eo(Ar) + A(iV)Ar^ ; e{N- At) = e^^N) + B{N)At'^ . (2.13) 

Acceptable fits were found for all datasets except for the E data at = 6, where the fit 
was restricted to Ar < 0.1. We find good evidence that the systematic error is indeed 
O(Ar^). Moreover, the coefficients A{N) and B{N) are themselves adequately fitted by the 
forms A{N) = aN"^, B{N) = bN"^ as displayed in Fig. |^. This behaviour is expected from 
inspection of the hybrid equations of motion [jT8| p9| . 

Since fluctuations in the system are suppressed by powers of 1/A^, we see that for 
= 0.25 systematic errors are likely to be dwarfed by statistical ones for reasonable sim- 
ulation runs. In the work presented in the rest of this paper we took Ar = 0.05 just to 
be conservative and safe. In fact, for very small or vanishing bare mass m near the critical 
point where the order parameter is very small, we also did production runs with Ar = 0.025 
to check explicitly that the systematic errors were smaller than the statistical errors. 

The bulk of the results of this paper were generated on a 16^ lattice with A^ = 0.25, 
using masses m ranging from 0.05 down to 0.0025. We have focussed our attention on the 
order parameter S, which differs from zero in the chiral limit only when chiral symmetry 

2 

is spontaneously broken, and is simply related to the chiral condensate: S = 
results for S as a function of bare mass m and inverse coupling are presented in Tab. 
|. Because the kinetic matrix M contains numerically large diagonal terms, this model is 
relatively cheap to simulate compared to, say, non-compact QED. This has enabled us to 
accumulate a large dataset, reflected in the relatively small statistical errors in Tab. |. 
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Typical runs which resulted in each of the entries in that table were between 15,000 and 
30,000 T units long. Such statistics are an order of magnitude better than present day state- 
of-the-art lattice QCD simulations. The statistical error bars in the table were obtained by 
binning methods which are particularly reliable when applied to such substantial data sets. 

In fact, because the interactions appear on the diagonal, M is sufficiently well-conditioned 
in the broken phase to permit simulations directly in the chiral limit m = 0. Of course, 
since spontaneous symmetry breaking cannot occur in a finite volume V, S is strictly zero 
in these simulations; defining 

* = ^E*(^). (2-14) 
then the next best thing to measure is 

|$| = {J.^itr$t|.y (2.15) 
This quantity is numerically very close to E extrapolated to the chiral limit, but exceeds it 



in a finite volume. We will discuss the finite volume correction in more depth in Sec. |V. 
Our results for |$| are given in Tab. |I[ To monitor finite volume effects directly we did a 
small number of runs on 12^ and 20'^ lattices, and tabulate the results in Tabs. |I| and |lll 



Finally in this section we give details of the model's gap equation. In the large- limit, 
for sufficiently strong coupling g^, the scalar auxiliary a develops a spontaneous vacuum 
expectation value E even in the chiral limit m ^ 0. To leading order in 1/A^ the relation 
between E, the bare mass m and is given by the lattice tadpole or gap equation 

1 m + E d'^k 1 

—r = l6 



c/2 E 7-^/2 (27r)4 sin^ k^+{m + E)2 

171 -\- / \ rv 

= ^rL ^«e^P-«(2 + (^ + S)')^'(f), (2.16) 

with Jq a modified Bessel function. In Fig. ^ we show predictions for E as a function of g"^ 
for m values 0.0, and 0.01. The m = line shows a continuous transition at l/g1 — 0.62 
between a symmetric phase at weak coupling and a broken phase at strong coupling - the 
curve having the {g"^ — g1)^ shape characteristic of mean field theory. For m > 0, E 7^ 
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for all values of g'^. Also shown are simulation results from a 6^ lattice with m = 0.01 for 
various A^, showing: 

(i) that increasing fluctuations, expected to be oc 1/iV, cause a suppression of S. 

(ii) that the simulated values of S exceed the gap equation predictions for > 3: this is 
in marked contrast to the case in three dimensions, where the gap equation gives an 
upper bound for all A^ [Q. Perhaps the lack of numerical accuracy of the gap equation 
is because the 1/A^ expansion is not renormalisable in four dimensions. 

III. NUMERICAL FITS TO THE EQUATION OF STATE 

In this section we report on fits to various trial forms for the equation of state m = 
/(S, l/g^) using data for the order parameter tabulated in Tab. |I[ We used the numerical 
package minuit to perform least squares fits, and in each case quote a standard error for the 
fitted parameters, and the P^i' degree of freedom, which gives a measure of the quality of 



First we list the functional forms we have examined, in increasing order of sophistication: 
I Mean Field (3 parameters): 



This is the simplest form, arising from a mean field treatment neglecting fluctuations 
in the order parameter field. 

II Power Law (5 parameters): 



This is a more general version of I, derivable by renormalisation group arguments from 
the assumption that a fixed point exists at the transition corresponding to a diverging 



the fit. 




(3.1) 




(3.2) 
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ratio between cutoff and physical scales. A useful discussion is found in [0]. The index 
6 is the standard critical exponent, and p = 6 — where /3 is another standard 
exponent. 

Ill Ladder Approximation (4 parameters): 



This is a restricted form of 11 based on the relation 7 = 1^5 — 1//3 = 1, inspired by 
solution of the quenched gauged U(l) NJL model in ladder approximation p5| . 

The last three fits all assume that the leading behaviour is described by the exponents 
of mean field theory (ie. S = 3, p = 1), but with logarithmic scaling corrections. Since 
these corrections introduce a dependence on some ultraviolet scale, they imply that the 
renormalised theory near the fixed point can never be made independent of the cutoff, and 
hence that the continuum limit is described by a free field theory - this is the phenomenon 
of triviality. The standard scenario is derived using renormalised perturbation theory in the 
context of 0(n) scalar field theory 0. This gives rise to: 

IV Sigma Model Logarithmic Corrections (4 parameters): 



Here the powers of the logarithmic terms are derived on the assumption that the 
effective field theory at the fixed point is the 0(4) linear sigma model, which has the 
same pattern of symmetry breaking as our SU(2)(8)SU(2) model. Note that we include 
C, the UV scale of the logarithm, as a free parameter. Since the condensate data are 
measured in lattice units, we expect a plausible fit should have C of 0(1). 

V Modified Sigma Model Logarithmic Corrections (5 parameters): 




(3.3) 




(3.4) 




(3.5) 
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Here we allow the powers of the logarithms to vary. This phenomenological form has 
been used extensively in fits of the equation of state of both non-compact lattice QED 
[|g and the U(l) NJL model |16[. Note that we have set the UV scale in the 



logarithm to 1; fits in which both p and C are allowed to vary are inherently unstable, 
since ln(C/S) = const. + ln^/^°'"(l/S), so that as S — the functional forms become 
so similar that the covariance matrix becomes singular. 

VI Large- Logarithmic Corrections (4 parameters): 

m = A(\-\]^ + BT.Hn (^] . (3.6) 



,9^ 9c/ 

This form is predicted for the NJL model in the large- limit P], and was used to fit 
data from the four dimensional Z2 NJL model in . There are important qualitative 
differences between form VI and the previous forms IV, V, due to the different role 
of the logarithmic term - i.e. the "effective" value of the exponent 6 measured at the 
critical coupling is larger than the mean field value 3 for IV and V, but smaller than 
3 for fit VI. 

We also tried to fit slightly modified forms of IV - VI, eg, allowing separate scales in 
the logarithms, but the results were either unstable, or not sufficiently distinct from the six 
forms presented here to be worth reporting. 

As we shall see, the main issue to arise when assessing the various forms ( p.l[ - |3.6|) is 
how much of the data to include in the fit. We expect any fit only to be successful in some 
finite scaling region around the transition. To see this, consider fits using all 84 data points, 
shown in Tab. |^ No fit gives acceptable results as judged by the large x^; although the 
two 5 parameter forms II and V are both noticeably better. To proceed, we must make 
some assumptions about the size and shape of the scaling window in the {l/g'^,m) plane, 
and truncate the dataset accordingly. In Tab. [V| we exclude extremal inverse coupling values 
and fit to G [0.52, 0.55], which still includes over half the data in the set. The resulting 
fits fall into two camps. Forms I, III and VI all yield similar large ^1 9l-i ^^^d are all 
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in effect forced very close to the mean field form 5 = 3, p = 1; in the case of VI by having 
the UV scale so large that the logarithmic term is effectively constant. These fits all fail 
because the data in this window does not favour p = 1. Forms II, V and to a lesser extent 
IV all allow the effective p to increase from 1, e.g. for form V peff = 1 + gi/ln(^) > 1. 
This enables better fits, which are also characterised by a smaller value for l/g^] to see why 
consider a Fisher plot of the data, in which is plotted against m/S. The Fisher plot is 
designed so that the mean field equation of state yields curves of constant l/g"^ which 
are straight lines of uniform slope, intercepting the vertical axis as m — in the broken 
phase, the horizontal axis as m — > in the symmetric phase, and passing through the origin 
as m ^ for g'^ = g1- Any departures from mean field behaviour show up as curvature in 
and variation in spacing between lines of uniformly spaced values of l/g"^. 

In Fig. ^we plot lines of constant 1/ g'^ using III, showing that the lines are straight and 
make a poor job of passing through the data points. The plots from I and IV are almost 
indistinguishable. In Fig.^ we show the same plot for II, and in Fig. |^ the plot for V. These 
last two fits are relatively successful because they accommodate lines of constant 1/ g'^ whose 
curvature changes sign according to which phase one is in. 

Outside the fitted window form II copes slightly better in the broken phase. It is worth 
noting that forms II and V gave fits of similar quality when applied both to simulation 
results of both the U(l) NJL model |jl6| and non-compact QED [|14|; clearly the logarithmic 



correction of ( |3.5|) can equally well be modelled by an effective 5 > 3. We have checked that 
fits II and V were stable under exclusion oil/g'^ = 0.52,0.55, and also under exclusion of 
the low mass points m = 0.005, 0.0025. 

A study of Figs. ||,|^ however, indicates that fits II and V are least satisfactory for the low 
mass points closest to the origin. This has led us to explore an alternative scaling window, 
in which all l/g"^ values are kept, but high mass points are excluded. Results from just 
keeping the 29 points corresponding to m = 0.0025, 0.005 are given in Tab. |VI|, and those 
for the 40 points obtained by also including m = 0.01 in Tab. IVTTl The picture changes 



dramatically: the values from fits III and VI are now much smaller and comparable to 

15 



those from fits II and V, though for m < 0.005 (Tab. |VlD the two 5 parameter fits still yield 
the lowest values. Once m = 0.01 points are included, fits II, III, V and VI are virtually 
indistinguishable in quality. 

What has happened is that the data from this window prefer an effective value of p 
much closer to 1, and an effective 6 now less than 3. This results, e.g. in the value of qi in 
fit V being very small, almost consistent with zero. The most impressive effect is that the 
logarithmic correction of form VI is now acting in the correct way to fit the data, moreover 
with a reasonable value for the UV scale C: in fit V this manifests itself in a relatively large 
negative value for the exponent q2. The form IV, which assigns a positive value to q2, is 
unstable. 

In Figs. P and ^ we show Fisher plots for the forms II and VI respectively, using the fit 
parameters of Tab. |VI1| , showing that the fits are indeed successful over a wide range of 1/(7^ 
at the expense of missing the higher mass points. The even spacing of the lines of constant 
l/g"^ in the broken phase is responsible for forcing the effective value of p to one. 

The values of per degree of freedom for the four successful fits are slightly too high for 
us to claim that they are the last word on understanding the model's equation of state. Even 
though the fits reproduce the spacing of the data in l/g'^, close inspection of Figs. ^ and |^ 
suggests that the variation of the data with m in the broken phase is not so well-described 
- clearly more low mass data will be needed to settle the issue, which in turn will require a 
quantitative understanding of finite volume corrections. We can claim, however, that there 
appears to be a dramatic switch in the preferred form of the equation of state according to 
whether the scaling window is chosen long and thin or short and fat in the {l/g'^,m) plane. 
In the next section we will consider data taken directly in the chiral limit, and find indirect 
evidence to support the short fat option. 
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IV. FINITE VOLUME CORRECTIONS IN THE CHIRAL LIMIT 



In the previous section we made no attempt to correct for finite volume effects in the 
measurements of S. From test runs for m = 0.005 dX 1/ = 0.50 on different-sized lattices, 
shown in Tab. |TT1|, finite volume effects are clearly present on a 12^ system, but the difference 
between 16^ and 20^, though statistically significant, is less than 2%. We hope, therefore, 
that finite volume corrections would have no impact on the qualitative conclusions of the 
previous section. 

Finite volume effects are numerically important, however, when considering measure- 
ments of the quantity I'I'I made with m = 0, defined in (|2.15|) and shown in Tab. ||. Naively, 
we expect that in a finite system |$| differs from the true order parameter S extrapolated 
to the chiral limit, because in the absence of a symmetry breaking term it is impossible to 
disentangle fluctuations of the order parameter field from those of the Goldstone modes; 
including the effects of the latter will give |$| > Sq, where Sq denotes the value of S in the 
chiral limit. Only in the thermodynamic limit V — > oo do we expect the Goldstone modes to 
average to zero and the two quantities to coincide. This phenomenon also occurs in numer- 



ical studies of scalar field theory, and has been analysed for the 0(4) sigma model in [jTT 



Formally, |$| corresponds to the minimum of a "constraint effective potential" in a finite 
volume V\ for the 0(4) sigma model this can be calculated using renormalised perturba- 
tion theory. There is indeed a correction due to the Goldstone modes which is intrinsically 
0{y^^). Numerically, however, a far more significant correction to the calculation arises 
from the need to renormalise the model in order to relate lattice parameters to measured 



quantities. To one-loop order this results in a relation [IT 

l-^-l = So (^1 + B{mRL)^^ + 0{gl)^ , (4.1) 

where m^j, qr are the renormalised scalar mass and coupling strength, L is the linear size of 
the system, and Z is the wavefunction renormalisation constant, defined as the coefficient 
of in the unrenormalised pion propagator in the limit — 0. The factor i? is a 
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slowly varying function of m/jL which accounts for the difference between one-loop integrals 
evaluated in finite and infinite volumes (for further background see [0). A fit to ( [4.1D for 
magnetisation data in the sigma model is given in Fig. 2 of [ pA]] . 

For the NJL model there is no reason for renormalised perturbation theory to apply; 
however since the effective theory in the broken phase is similar to the sigma model, in the 
sense that it has the same light scalar degrees of freedom, we expect a result similar to ( |4.1|) 
to hold. We will therefore attempt to account for the finite volume correction A = |$| — Sq 
with the formula 

A = (4.2) 

where the B factor has been set constant, and higher order effects are ignored. 

The first prediction of ( [4.2|) is that A oc L^^. To test this, we can compare |$| data taken 
at 1/(7^ = 0.48 from three different volumes against our best estimate for Sq. From the data 
of Tab. we see that |$| decreases with L as predicted. To evaluate Sq, we use the two 
most successful fits to the equation of state from the previous section, each one characteristic 
of the two different hypotheses about the shape of the scaling window. The narrow scaling 
window 0.52 < l/g"^ < 0.55 was best fitted by form V ( ^.51 ) using the parameters of Tab. |V| 
(the corresponding Fisher plot is shown in Fig. in the chiral limit at = 0.48 this 
form predicts Sq = 0.3341 (we implicitly assume here that the equation of state fits apply to 
the thermodynamic limit). This is already equal to the value of I'I'I on a 12^ lattice, and lies 
above |$| on larger systems, making it difficult to see how ( ^4.2| ) applies. The second estimate 
comes from assuming a wide scaling window in l/g"^ but excluding higher mass data. This 
was best described using form VI ( |3.6| ) using the parameters of Tab. [VII| (Fig. H); the chiral 



limit prediction now is Sq = 0.3030. In Fig. |10| we plot (y|)"A vs. 1/L for n = 1, 2 and 3 
and L = 12, 16 and 20, assuming Sq = 0.3030. We find that L^A is approximately constant, 
which supports the hypothesis ( [4.2| ). 

Next we examine the dependence of A on 1/ (7^, and hence So, using data from the broken 
phase of the model taken from a constant volume of 16^. The data of Tab || were used in 



conjunction with the same two fits for Sq to produce the two sets of values for A given 
in Tab. |V11J|. The quoted errors are the statistical errors in the measurement of |$|, and 



take no account of errors in extrapolating S to the chiral limit. The values Ay obtained 
using form V actually change sign over the range of l/g"^ explored, and clearly cannot be 
fitted by any relation of the form ( [4.21 ). The values Ayj obtained using form VI are plotted 



against 1/ g'^ in Fig. |TT|. The most interesting trend is that Ayj is approximately constant 
for 0.45 < l/g"^ < 0.52, whereas So itself falls from 0.4115 to 0.1332 over the same range. 
How do we reconcile this observation with (|4.2| ), which apparently predicts A oc I/Sq? 



The answer lies in the field dependence of the wavefunction renormalisation constant Z, 
which differs markedly between the ferromagnetic 0(4) sigma model and the fermionic NJL 
model 1^ [|1^. In the sigma model, Z is perturbatively close to 1 (as borne out by the 
data of |jll|]), and hence A oc I/Sq. In the NJL model, on the other hand, the large- Nf 
approximation predicts Z oc l/ln(l/So) [§, and hence 

A oc \,y (4.3) 

We can see this using another argument. In both sigma and NJL models, a Ward identity 
plus the assumption that both conserved current and transverse field couple principally to 
the Goldstone mode (in particle physics this mode is the pion, and the assumption is PCAC), 
leads to the relation p6[ 

Zf^ = (4.4) 

where f-j^ is the coupling to the axial current, or in the context of hadronic physics the pion 
decay constant. We thus have 

Now, in the large- A^/ approximation (see Appendix for an explicit calculation), 

./AjSo 1 / 1 \ 
/. = ^^^lnM^)> (4-6) 
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yielding once again the relation ( |4.3| ). 

Note that both scenarios predict triviality of the resulting theory. In the NJL model, in 
the large- A^/ treatment Sq is a physical scale related to the renormalised fermion and scalar 
masses. The pion decay constant in physical units is thus /tt/Sq, which diverges as ln2 (I/Sq) 
in the continuum limit Sq 0. In the sigma model, /tt/Sq = which is constant in the 
continuum limit; however the physical scale is now the scalar mass m„ oc ^/gflIlQ. Hence in 
physical units 

■^'^ oc ^ ocln5 (^], (4.7) 



and again diverges in the continuum limit 



We tested the hypothesis ( [4.3| ) with a two parameter fit of the data of Tab. [VIII| to the 
form 

A = ^ , (4.8) 



.So, 

With all points included we find A = 1.86(13), B = 1.10(11) with xVdof of 4.3; however if 
the l/g"^ = 0.53 point is left out the fit is much more satisfactory: A = 1.40(10), B = 0.82(6) 
with x^/dof of 0.81. This fit is plotted in Fig. [Tl|. In either case it is reassuring that the 
UV scale of the logarithm is close to unity. It is also interesting that the fit accounts for 
the slight curvature in the plot, though this aspect is not tightly constrained by the current 
error bars. The fitted curve also rises sharply to pass close to the excluded point. 

To conclude, we can account for our finite volume effects in the chiral limit under three 
assumptions: 

(i) that the equation of state fits to the S data from the 16^ lattice in the previous section, 
when extrapolated to the chiral limit, accurately describe the thermodynamic limit. 



(ii) that the most appropriate form to extrapolate the equation of state is form VI ( p^ 



using the parameters of Tab. |VII . 

(iii) that the finite volume correction has the form (f4.3|) 
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The latter two assumptions are consistent with logarithmic corrections of the form ad- 
vocated in , which differ qualitatively from those used to describe triviality in scalar field 
theory. 



V. DISCUSSION 

It would appear that for the lattice sizes and bare masses used in this study, which are 
close to state-of-the-art (though recall that dual-site four-fermi models are relatively cheap to 
simulate), it is still difficult to distinguish different patterns of triviality, or indeed distinguish 
a trivial from a non-trivial fixed point, purely on the basis of fits to model equations of state; 
different assumptions about the size and shape of the scaling region result in different best 
fits. Only once the analysis of Sec. IV is applied to the data taken at zero bare mass does 
the form ( p.6[ ), associated with composite scalar degrees of freedom, appear preferred. Our 
principal conclusion is thus that the form ( |1.4| ), ( |3.6| ) is qualitatively correct beyond leading 
order in the 1/iV/ expansion. 

The sensitivity of the preferred form of the equation of state to the shape of the assumed 
scaling window has implications for similar studies in related fermionic models, noteably the 
U(l) NJL model studied using the "link" formulation [Q, and non-compact QED4 fl^; the 



two recent studies cited attempt fits based on the power law form (|3.2|) and the modified 
sigma model form (|3.5| ), and find the two forms difficult to distinguish on the basis of 
alone (in other words, a value of 5 > 3 is difficult to distinguish from a logarithmic correction 



with q2 > 0). In Tab. |^ the values of qi and q2 are shown for the two scaling windows 



used in this study, together with the quoted fits for the two other models from [T^ and 



as well as theoretical values for the 0(4) and 0(2) sigma models [^] [0, corresponding 
respectively to the broken SU(2) C?>SU(2) symmetry in this study, and the broken U(l) chiral 
symmetry of the models of and Using the data of one can examine the stability 



of the fitted qi and q2 with respect to exclusion of higher mass data; the quoted value for qi 
appears quite stable, that for q2 less so. 
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An outstanding issue is whether the disparity in the fitted values of the q exponents in 
Tab. |I^ means that the different lattice models lie in different universality classes. The 
analysis presented in this paper suggests that current lattice simulations are unable to de- 
cide. Unfortunately for the two models considered in [|T6| there is no realisable way of 
simulating directly in the chiral limit, and hence the analysis of Sec. IV is unavailable. A 
useful direction to explore may be a simulation of the U(l) NJL model with the dual site 
formulation. 
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APPENDIX: l/Np CALCULATION OF 

In this appendix we calculate the pion decay contant /^r, defined by the PC AC relation 
(0|yl^|7r°(g)) = {Q\^^^^^Jl^\^-{q)) = rq^f^^ (Al) 
in the large- A^/ approximation of the NJL model. Diagrammatically the matrix element to 



analyse is shown in Fig. 12. Note that the external pion leg yields a normalisation factor 



Z2, where Z is defined by the infrared limit of the pion propagator Dt^: 



Z 

q^O q 



limD^iq) = -. (A2) 



For d < A Z is finite. In the NJL model to leading order in 1/A^/ it is given by 



1 (An) 2 
?2r(2-f)' 



Z = 3;^TT7T7^S^-'^- (A3) 



The matrix element is then given by 
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9 



tr 



1 



1 



-7/^75^ 



-75^ 



(A4) 



where details of Feynman rules and evaluation of are also given in [Q. For d < 4 the 



integral can be evaluated to yield 

{0\A''^\7r"{q))=gZ-2jN'fX4Eiq^x 



r(2-f^ /•! 



dx- 



(47r)2 ^0 (S2 + 



2\2- 



(A5) 



We now take the chiral limit ^ (f-, valid for an on-shell pion; combining this step with 



relations (|A1|) and (|A3|) we find 



r(2 



\ (47r)f 



(A6) 



This agrees with the result found for = 3 in ||27|| . 

Finally we examine the limit d ^ 4. Defining e = 4 — (i, we have 



'Nf^ 1 1, S2 



IE 



TT 



' In — + 0(e). 

e 2 4:71 2 ^ ' 



(A7) 



In = 4 the model is no longer renormalisable, and we must introduce an explicit UV cutoff, 
which can be set to one by appropriate choice of units. On the assumption that the 1/e 
term is then replaced by the logarithm of the cutoff, we are left with the scaling prediction 



ln5 -. 

vr L 



(A8) 
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TABLES 



TABLE I. Order parameter data S from a 16^ lattice. 



1 / 2 


m = 0.05 


m = 0.04 


m = 0.03 


m = 0.02 


m = 0.01 


m = 0.005 


m = 0.0025 


0.45 


0.5503(3) 


0.5344(3) 


0.5149(3) 


0.4919(4) 


0.4611(9) 


0.4411(9) 


0.4241(13) 


0.46 












0.4080(11) 


0.3873(18) 


0.47 












0.3753(10) 


0.3509(19) 


0.48 












0.3452(11) 


0.3213(15) 


0.49 












0.3091(12) 


0.2809(21) 


0.50 


0.4371(4) 


0.4152(5) 


0.3899(6) 


0.3570(7) 


0.3091(7) 


0.2746(8) 


0.2470(18) 


0.51 












0.2435(8) 


0.2107(19) 


0.52 


0.3970(3) 


0.3733(5) 


0.3443(6) 


0.3086(5) 


0.2520(5) 


0.2097(9) 


0.1774(20) 


0.525 


0.3877(3) 


0.3636(3) 


0.3337(3) 


0.2954(3) 


0.2398(4) 


0.1928(8) 


0.1583(17) 


0.53 


0.3780(3) 


0.3530(3) 


0.3228(3) 


0.2836(5) 


0.2257(7) 


0.1789(8) 


0.1425(19) 


O.OOO 




O.o43o(4j 


n Qi on/' o\ 
U.ol29(oj 


U.2l6\J[0) 


U.zi4z(^D j 


0.16o/(9j 


0.12oo(20j 


0.54 


0.3592(2) 


0.3329(2) 


0.3011(3) 


0.2611(7) 


0.2004(6) 


0.1506(11) 


0.1081(16) 


0.545 


0.3505(3) 


0.3243(3) 


0.2923(3) 


0.2504(3) 


0.1879(5) 


0.1355(10) 


0.0926(20) 


0.55 


0.3409(3) 


0.3152(3) 


0.2825(3) 


0.2405(7) 


0.1773(9) 


0.1246(13) 




0.60 


0.2641(3) 


0.2353(3) 


0.2012(6) 


0.1568(7) 


0.0941(5) 


0.0526(10) 




0.70 


0.1652(3) 


0.1396(3) 


0.1112(4) 


0.0784(4) 


0.0403(3) 


0.0207(4) 
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TABLE II. Measured values of |$| for m = from various lattice sizes. 



1/^2 12^ 16^ 20"^ 

0.45 - 0.4301(11) 

0.46 - 0.3955(11) 

0.47 - 0.3587(8) 

0.48 0.3341(11) 0.3213(9) 0.3164(8) 

0.49 - 0.2830(10) 

0.50 - 0.2424(10) 

0.51 - 0.2013(14) 

0.52 - 0.1574(14) 

0.53 - 0.1149(18) 



TABLE III. Order parameter S for m = 0.005, l/g^ = 0.50 from various lattice sizes. 

12^ 16^ 20^ 

0.2643(14) 0.2746(8) 0.2786(7) 
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TABLE IV. Equation of state fits using all data (84 points). 

A B \/gl 5 p (^^ C xVdof 

1^2 

I 1.8608(17) 0.8443(10) 0.52812(6) - - - - 559.2 

II 2.3757(74) 0.7636(27) 0.53941(25) 2.6901(68) 1.1978(21) - - 73.0 

III 1.7974(18) 0.5810(10) 0.55759(27) 2.2237(40) - - - 183.4 

IV 2.884(13) 1.497(18) 0.52006(9) - - - 2.850(59) 385.7 

V 1.9669(26) 0.9590(22) 0.53311(14) - - nn^IoS! ^ 39.7 

— 0.0592(52) 

VI 1.7877(17) 0.6577(33) 0.54635(12) - - - 2.004(12) 182.6 



TABLE V. Equation of state fits using data for 0.52 <l/g'^ < 0.55 (48 points). 

A B 1/gl 5 p C xVdof 

[02 

I 1.8152(53) 0.9616(13) 0.53309(7) - - - - 45.0 

II 3.845(71) 1.309(20) 0.52224(44) 3.50(2) 1.59(1) - - 3.1 

III 1.8160(53) 0.9943(89) 0.53228(22) 3.044(11) - - - 45.0 

IV 3.285(58) 2.64(11) 0.5274(2) - - - 7.13(82) 32.0 

V 2.160(10) 0.8499(53) 0.52562(32) - - 372(17) " ^'^ 

VI 1.8146(53) 0.0258(5) 0.53359(8) - - - 0.70(45) xlO^^ 45.0 
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TABLE VI. Equation of state fits using data for m < 0.005 (29 points). 





A 


B 


V9l 


5 p 


{ 

U2 


C 


xVdof 


I 


1.136(13) 


0.4907(73) 


0.51910(33) 








41.0 


II 


1.099(35) 


0.312(11) 


0.53781(83) 


2.226(26) 0.921(10) 






1.5 


III 


1.367(19) 


0.3992(64) 


0.53779(83) 


2.351(19) 






3.4 


IV 
V 

VI 


1.205(33) 
1.354(18) 


0.633(18) 
0.421(16) 


0.53157(54) 
0.53334(52) 


no fit found 


-0.118(25) 
-0.707(34) 


2.01(8) 


2.4 
2.9 


TABLE VII. Equation of state fits using data for m < 0.01 (40 points). 




A 


B 




S p 




C 


xVdof 


I 


1.403(6) 


0.6349(43) 


0.52297(18) 








72.4 


II 


1.382(17) 


0.3998(68) 


0.54149(61) 


2.267(18) 0.9731(53) 






3.9 


III 


1.465(7) 


0.4316(35) 


0.54019(83) 


2.334(12) 






4.5 


IV 








no fit found 








V 
VI 


1.476(13) 
1.459(7) 


0.7837(73) 
0.474(10) 


0.53391(32) 
0.53546(31) 




0.016(11) 
-0.553(18) 


1.95(5) 


4.6 
4.0 
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TABLE VIII. Values of the quantity A for the two fits to the thermodynamic and chiral hmit 
discussed in the text. 



1/5^ 






^V 


^VI 




error 


0.45 






-0.0326 


0.0186 




0.0011 


0.46 






-0.0253 


0.0196 




0.0011 


0.47 






-0.0195 


0.0189 




0.0008 


0.48 






-0.0128 


0.0183 




0.0009 


0.49 






-0.0045 


0.0181 




0.0010 


0.50 






0.0058 


0.0175 




0.0010 


0.51 






0.0233 


0.0195 




0.0014 


0.52 






0.0568 


0.0242 




0.0014 


0.53 






0.1149 


0.0436 




0.0018 




TABLE IX. 


Values of the q exponents for the various models discussed in Sec. 


V. 




0(4) 


0(2) 


SU(2) NJL 


SU(2) NJL 


U(l) NJL 


QED4 








(Tab. 0) 


(Tab. |VI| 


(Ref. @) 


Ref. 0) 


Qi 


0.5 


0.4 


0.786(19) 


0.016(11) 


0.36(11) 


0.485(7) 


<12 


1 


1 


0.372(17) 


-0.553(18) 


0.84(1) 


0.324(15) 
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FIG. 1. The condensate S plotted against timestep At for N = 0.25 (o), iV = 1 (A), N = 3 
(□) and iV = 6 (o) 
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FIG. 2. Same as Fig. |T] but for the energy density e 



32 



4.0 



3.0 



2.0 



OA 
□ B 



1.0 



0.0 -Q-= = - 



-1.0 



1 



3 4 
N 



\ 

6 



FIG. 3. Behaviour of the coefficients A{N) and B{N) 
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FIG. 4. Gap equation predictions for S for m = (sohd hne) and m = 0.01 (dashed hne), 
together with simulation data (see caption to Fig. |l|) 
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FIG. 8. Fisher plot for fitting form II, m < 0.01. 
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FIG. 9. Fisher plot for fitting form VI, m < 0.01. 
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FIG. 10. (L/12)"A vs. 1/L for Xjg^ = 0.48. 
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